#include "bcpp_read_seq.h"
#include <boost/regex.hpp>
#include <boost/tuple/tuple.hpp>
#include <boost/lexical_cast.hpp>
#include <boost/foreach.hpp>
#include <boost/format.hpp>
#include <iostream>
#include <algorithm>
#include <stdexcept>
#include "blast_results.h"
#include "ultoa.h"
#include "print_error.h"
#include "bcpp_parsefasta.h"
#include "bcpp_parseclustal.h"
#define _MY_STL_MACROS
#include "stlhelper.h"
#include "cdna_prot.h"
#include "print_sequence.h"

using std::string;
using std::vector;

namespace BioCPP
{
typedef std::vector<unsigned>::const_iterator cu_iter;
//________________________________________________________________________________________

//	t_collect_sequences_from_strm

//________________________________________________________________________________________
struct t_collect_sequences_from_strm
{
	std::string							seqid_regex_str;
	t_progress_indicator&				dots;
	t_collect_sequences_report_errors& 	error_report;
	cdna_prot::t_validate_prot			valid_seq_char;
	unsigned							acc_beg_line;	// accession begin line for error

	// data to be retrieved
	t_collect_sequences&				seq_db;

	t_collect_sequences_from_strm(  const std::string&				seqid_regex_str_,
									t_progress_indicator&			dots_,
									t_collect_sequences_report_errors& error_report_,
									const std::string&				valid_sequence_chars,
									bool							to_uppercase,
                                    bool                            letters_are_valid,
                                    t_collect_sequences&			seq_db_)
			:
			seqid_regex_str(seqid_regex_str_),
			dots(dots_),
			error_report(error_report_),
			valid_seq_char(valid_sequence_chars, to_uppercase, false, letters_are_valid),
			acc_beg_line(0),
			seq_db(seq_db_)
		 {	}



	// callback
	bool operator()(std::string& accession, std::string& sequence,
						unsigned acc_filepos,
						unsigned seq_filepos);
    // callback
    bool operator()(std::string& accession, std::string& sequence,
                        unsigned curr_line);
};


// 
// 
//  line counts for error messages
// 
bool t_collect_sequences_from_strm::operator()(std::string& accession, std::string& sequence, unsigned acc_beg_line)
{
	// print progress dots
	++dots;

	// remove invalid characters
	valid_seq_char.validate(sequence);

	// bad sequence!!!
	if (!sequence.length())
	{
		error_report.zero_length_sequence(acc_beg_line);
		return true;
	}
	if (!accession.length())
	{
		error_report.zero_length_accession(acc_beg_line);
		return true;
	}

	// match id
    static boost::regex seqid_regex(seqid_regex_str);
    boost::match_results<std::string::iterator> id_match;
    string sequence_id;

    // fail to match regex => error
    if (!boost::regex_search(accession.begin(), accession.end(), id_match, seqid_regex) &&
        id_match[1].length())
    {
        error_report.no_matching_seq_id(acc_beg_line, seqid_regex_str);
        return true;
    }

    // is this a new accession?
    bool new_seqid_added = false;
    unsigned name_index = seq_db.seq_names.index(id_match[1], new_seqid_added);
    if (new_seqid_added)
    {
        seq_db.seq.push_back(sequence);
        seq_db.acc.push_back(accession);
        seq_db.name_indices_to_seq_indices.push_back(name_index);
    }

    // this accession has been used
    else
    {
        // allow duplicates
        if (seq_db.allow_duplicates == eSEQ_DUPLICATES)
        {
            seq_db.seq.push_back(sequence);
            seq_db.acc.push_back(accession);
            seq_db.name_indices_to_seq_indices.push_back(name_index);
        }

        // append to sequence
        else if (seq_db.allow_duplicates == eSEQ_APPENDDUPLICATES)
            seq_db.seq[name_index] += sequence;

        // error
        else
            error_report.duplicated_seq_id(name_index, acc_beg_line);
    }

	return true;
}


bool t_collect_sequences_from_strm::operator()(std::string& accession, std::string& sequence,
									unsigned /*acc_filepos*/, unsigned /*seq_filepos*/)
{
	// print progress dots
	++dots;

	// get line counts for error messages
	unsigned curr_acc_beg_line = acc_beg_line;
	acc_beg_line += std::count(sequence.begin(), sequence.end(), '\n') + 1;

	// remove invalid characters
	valid_seq_char.validate(sequence);

	// bad sequence!!!
	if (!sequence.length())
	{
		error_report.zero_length_sequence(curr_acc_beg_line);
		return true;
	}
	if (!accession.length())
	{
		error_report.zero_length_accession(curr_acc_beg_line);
		return true;
	}

	// match id
    static boost::regex seqid_regex(seqid_regex_str);
    boost::match_results<std::string::iterator> id_match;
    string sequence_id;

    // fail to match regex => error
    if (!boost::regex_search(accession.begin(), accession.end(), id_match, seqid_regex) &&
        id_match[1].length())
    {
        error_report.no_matching_seq_id(curr_acc_beg_line, seqid_regex_str);
        return true;
    }

    // is this a new accession?
    bool new_seqid_added = false;
    unsigned name_index = seq_db.seq_names.index(id_match[1], new_seqid_added);
    if (new_seqid_added)
    {
        seq_db.seq.push_back(sequence);
        seq_db.acc.push_back(accession);
        seq_db.name_indices_to_seq_indices.push_back(name_index);
    }

    // this accession has been used
    else
    {
        // allow duplicates
        if (seq_db.allow_duplicates == eSEQ_DUPLICATES)
        {
            seq_db.seq.push_back(sequence);
            seq_db.acc.push_back(accession);
            seq_db.name_indices_to_seq_indices.push_back(name_index);
        }

        // append to sequence
        else if (seq_db.allow_duplicates == eSEQ_APPENDDUPLICATES)
            seq_db.seq[name_index] += sequence;

        // error
        else
            error_report.duplicated_seq_id(name_index, curr_acc_beg_line);
    }

	return true;
}


//________________________________________________________________________________________

//	t_collect_sequences_report_errors_to_strm

//________________________________________________________________________________________
void t_collect_sequences_report_errors_to_strm::zero_length_sequence(unsigned line)
{
	CERR << "Error:\n\tEmpty sequence at line " << line << "\n";
}

void t_collect_sequences_report_errors_to_strm::zero_length_accession(unsigned line)
{
	CERR << "Error:\n\tEmpty accession at line " << line << "\n";
}

void t_collect_sequences_report_errors_to_strm::no_matching_seq_id(unsigned line, const string& regex_str)
{
	CERR << "Error:\n\tNo matching sequence id for the sequence at line " << line <<
					" with the regular expression: \"" << regex_str << "\"\n";
}

void t_collect_sequences_report_errors_to_strm::duplicated_seq_id(unsigned previous_index, unsigned line)
{
	CERR << "Error:\n\tThe sequence id at line " << line << " has already been used by sequence #" <<
					previous_index + 1 << ", and will be ignored\n";
}


//________________________________________________________________________________________

//	t_collect_sequences_collect_errors

//________________________________________________________________________________________
void t_collect_sequences_collect_errors::report(std::ostream& os, bool detail) const
{
	if (empty_sequence_line_num.size())
	{
		std_print_error(os, ultoa(empty_sequence_line_num.size()) +
							"\tempty sequences were found (Accessions without sequence)");
		if (detail)
		{
			os << "\ton the following lines:\n";
			BOOST_FOREACH(unsigned i, empty_sequence_line_num) 
				os << "\t" << i << "\n";
		}
	}
	if (empty_accession_line_num.size())
	{
		std_print_error(os, ultoa(empty_sequence_line_num.size()) +
							"\tmissing accessions were found "
							"(Sequences without accession lines)");
		if (detail)
		{
			os << "\ton the following lines:\n";
			BOOST_FOREACH(unsigned i, empty_accession_line_num) 
				os << "\t" << i << "\n";
		}
	}
	if (no_matching_seq_id_line_num.size())
	{
		std_print_error(os, ultoa(no_matching_seq_id_line_num.size()) +
							"\tAccessions did not have identifiers matching the supplied "
							"regular expression \"" + regex_str + "\" and were discarded\n");
		if (detail)
		{
			os << "\ton the following lines:\n";
			BOOST_FOREACH(unsigned i, no_matching_seq_id_line_num) 
				os << "\t" << i << "\n";
		}
	}
	if (sequence_id_duplicated_line_num.size())
	{
		std_print_error(os, ultoa(sequence_id_duplicated_line_num.size()) +
							"\tDuplicate sequence identifiers were found previously "
							"and the sequencs were discarded.\n");
		if (detail)
		{
			os << "\tThe indices of the originals and the line numbers of the duplicates:\n";
			for (unsigned i = 0; i < sequence_id_duplicated_line_num.size(); ++i)
			{
				os	<< "\t"
					<< sequence_id_duplicated_prev_index[i] << "\t"
					<< sequence_id_duplicated_line_num[i] 	<< "\n";
			}
		}
	}
}
void t_collect_sequences_collect_errors::zero_length_sequence(unsigned line)
{
	empty_sequence_line_num.push_back(line);
}

void t_collect_sequences_collect_errors::zero_length_accession(unsigned line)
{
	empty_accession_line_num.push_back(line);
}

void t_collect_sequences_collect_errors::no_matching_seq_id(unsigned line, const string& regex_str_)
{
	no_matching_seq_id_line_num.push_back(line);
	regex_str = regex_str_;
}

void t_collect_sequences_collect_errors::duplicated_seq_id(unsigned previous_index, unsigned line)
{
	sequence_id_duplicated_prev_index.push_back(previous_index);
	sequence_id_duplicated_line_num  .push_back(line);
}





//________________________________________________________________________________________

//	t_collect_sequences

//________________________________________________________________________________________
void t_collect_sequences::init_from_stream(	std::istream&					istrm,
											std::string						regex_str,
											t_progress_indicator&			dots,
											t_collect_sequences_report_errors& errors,
											const std::string&				letters,
											bool							to_uppercase,
                                            bool                            letters_are_valid)
{
	boost::regex brackets("\\(.*\\)");
	if (!boost::regex_search(regex_str, brackets))
		regex_str = "(" + regex_str + ")";

	t_collect_sequences_from_strm collect_seq
									( regex_str, dots, errors, letters, to_uppercase, letters_are_valid, *this);
    try
    {
        BioCPP::Process_FA_Sequences_tellg(istrm, collect_seq);
    } 
    catch (clustal_format_error& error)
    {
        set_allow_duplicates(eSEQ_APPENDDUPLICATES);
        BioCPP::Process_Clustal_Sequences(istrm, error.preamble, error.matching_format_line,
                                          error.line_cnt, collect_seq);
    }
}



bool t_collect_sequences::sequence(const string& name, string&  seq) const
{
	unsigned name_index;
	if(!get_index_of_name(name, name_index))
		return false;
	seq = sequence(name_index);
	return true;
}
void t_collect_sequences::sequence(const std::string& name, std::vector<string>&  seq) const
{
	vector<unsigned> indices;
	get_indices_of_name(name, indices);
	vector<unsigned>::iterator i = indices.begin();
	vector<unsigned>::iterator end = indices.end();
	while (i != end)
	{
		seq.push_back(sequence(*i));
		++i;
	}
}
bool t_collect_sequences::sequence(const string& name, string&  seq, string& acc) const
{
	unsigned seq_index;
	if(!get_index_of_name(name, seq_index))
		return false;
	seq = sequence(seq_index);
	acc = accession(seq_index);
	return true;
}

void t_collect_sequences::sequence(const std::string& name, std::vector<string>&  seq,
										std::vector<string>& acc) const
{
	vector<unsigned> indices;
	get_indices_of_name(name, indices);
	vector<unsigned>::iterator i = indices.begin();
	vector<unsigned>::iterator end = indices.end();
	while (i != end)
	{
		seq.push_back(sequence(*i));
		acc.push_back(accession(*i));
		++i;
	}
}
bool t_collect_sequences::get_index_of_name(const std::string& name, unsigned& name_index) const
{
	// name not found
	if (!seq_names.is_indexed(name, name_index))
		return false;

	// no duplicates, order of seq names is the order of names
	if (allow_duplicates != eSEQ_DUPLICATES)
		return true;

	// find first in name_indices_to_seq_indices looking for name_index
	std::vector<unsigned>::const_iterator i =
					std::find(name_indices_to_seq_indices.begin(), 
                              name_indices_to_seq_indices.end(), name_index);
	assert(i != name_indices_to_seq_indices.end());
	name_index = i - name_indices_to_seq_indices.begin();
	return true;
}

void t_collect_sequences::get_indices_of_name(const std::string& name,
												std::vector<unsigned>& indices) const
{
	unsigned name_index;

	// name not found
	if (!seq_names.is_indexed(name, name_index))
		return;


	// no duplicates, order of seq names is the order of names
	if (allow_duplicates !=  eSEQ_DUPLICATES)
	{
		indices.push_back(name_index);
		return;
	}

	// find all in name_indices_to_seq_indices which match
	std::vector<unsigned>::const_iterator i = std::find(name_indices_to_seq_indices.begin(),  
                                                        name_indices_to_seq_indices.end(), name_index);
	while (i != name_indices_to_seq_indices.end())
	{
		indices.push_back(i - name_indices_to_seq_indices.begin());
		i = std::find(++i, name_indices_to_seq_indices.end(), name_index);
	}
	return;
}

void t_collect_sequences::replace_names(const t_index_string&	new_seq_names)
{	
	assert(seq_names.size() ==new_seq_names.size());
	seq_names = new_seq_names;

	 for (unsigned i = 0; i < acc.size(); ++i)
		 acc[i] = name(i);
}
void t_collect_sequences::print_fasta(	std::ostream&					os,
										const std::vector<unsigned>&	valid_indices,
										unsigned						line_wrap) const

{
	// indices range
	cu_iter	index = valid_indices.begin();
	cu_iter	end = valid_indices.end();

	for (; index  != end; ++index)
			::print_fasta(	os,
							accession(*index),
							sequence(*index),
							line_wrap);
}

void t_collect_sequences::print_phylip(std::ostream&					os,
										const std::vector<unsigned>&	valid_indices,
										unsigned						line_wrap,
										bool							interleaved) const

{
	// print count and length
	os << valid_indices.size() << " " << sequence(valid_indices[0]).length() << "\n";

	print_multiple_aligned(os, valid_indices, line_wrap, interleaved, true, false);
}



void t_collect_sequences::print_clustal(std::ostream&					os,
										const std::vector<unsigned>&	valid_indices,
										unsigned						line_wrap) const

{
	os << "CLUSTAL W (1.8) multiple sequence alignment\n\n";

	print_multiple_aligned(os, valid_indices, line_wrap, true, false);
}


//
//	helper for phylip and clustal
//
namespace
{
struct t_invalid_aa_to_hyphens
{
    char table [256];

    // if stop_is_X, stop asterisks are replaced by X
    t_invalid_aa_to_hyphens()
    {
        std::fill_n(table, 256, '-');
        for (char c = 'a'; c < 'z'; ++c)
            table[static_cast<int>(c)] = c;
        for (char c = 'A'; c < 'Z'; ++c)
            table[static_cast<int>(c)] = c;
    }
    void validate(std::string& sequence)
    {
        std::string::iterator first	=  sequence.begin();
        std::string::iterator last	=  sequence.end();
    
        for (; first != last; ++first)
            *first = table[static_cast<int>(*first)];
    }
};
};
void t_collect_sequences::print_multiple_aligned(std::ostream&					os,
													const std::vector<unsigned>&	valid_indices,
													unsigned						line_wrap,
													bool							interleaved,
													bool							names_only_first_line,
                                                    bool                            invalid_to_underscores) const
{
	// indices range
	cu_iter	index = valid_indices.begin();
	cu_iter	end = valid_indices.end();


	// check all the same length
    unsigned initial_sequence = *index;
	unsigned seq_len = sequences()[initial_sequence].length();
	while (++index != end)
	{
		if (sequences()[*index].length() != seq_len)
        {
            string err_msg = (boost::format("sequence[%1%] = %2%, sequence[%3%] = %4%")
                                                        % initial_sequence
                                                        % seq_len 
                                                        % unsigned (*index) 
                                                        % sequences()[*index].length()).str();
            throw std::runtime_error(	"Sequences to be printed in the phylip format must "
                                        "be of the same length.\n" + err_msg + "\n");
        }                            
	}


	if (line_wrap == string::npos)
		line_wrap = seq_len;


	//
	// make sure all names same length padded with spaces. Should be 10 but who cares
	//
	// minimum padded length = 10
	std::size_t longest_name_len = 10;
	index = valid_indices.begin();
	while (index != end)
	{
		longest_name_len = std::max(longest_name_len, name(*index++).length());
	}
	// add two spaces between name and sequence
	longest_name_len+=2;

    t_invalid_aa_to_hyphens invalid_to_hyphens;

	//
	// print sequences
	//
	if (interleaved)
	{
		// only print names once if specified otherwise print_names always positive
		// even if line wrap once per letter
		int print_names = names_only_first_line ? 1 : seq_len;
		unsigned pos = 0;
		while (pos < seq_len)
		{
			for (index = valid_indices.begin(); index != end; ++index)
			{
				// print names
				if (print_names > 0)
				{
					os << name(*index);
					os << string(longest_name_len - name(*index).length(), ' ');
				}

                string seq_line = sequence(*index).substr(pos, line_wrap);
                if (invalid_to_underscores)
                    invalid_to_hyphens.validate(seq_line);

				// print sequence
				os << seq_line << "\n";
			}
			pos += line_wrap;
			print_names--;
			// extra line between sequences
			os << "\n";
		}
	}

	//
	// non-interleaved
	//
	else
	{
		for (index = valid_indices.begin(); index != end; ++index)
		{
			// print names
			os << name(*index);
			os << string(longest_name_len - name(*index).length(), ' ');

			// print remainder of sequence in non-interleaved format
			unsigned pos = 0;
			while (pos < seq_len)
			{
                string seq_line = sequence(*index).substr(pos, line_wrap);
                if (invalid_to_underscores)
                    invalid_to_hyphens.validate(seq_line);

                // print sequence
                os << seq_line << "\n";

				pos += line_wrap;

			}
		}
	}
}



void t_collect_sequences::print_sequences(	eSEQ_FORMAT 			format,
											std::ostream&			os,
											const vector<unsigned>& valid_indices,
											unsigned				line_wrap) const
{
	// no line wrap if 0
	if (!line_wrap)
		line_wrap = unsigned(-1);

	if (valid_indices.size() == 0)
		throw std::logic_error("No sequences to be printed.");


#ifdef DEBUG
	for (unsigned i = 0; i < valid_indices.size(); ++i)
		assert(valid_indices[i] < size());
#endif

	switch (format)
	{
		case eSF_FASTA:
			print_fasta(os,  valid_indices, line_wrap);
			break;
		case eSF_CLUSTAL:
			print_clustal(os,  valid_indices, line_wrap);
			break;
		case eSF_PHYLIP:
			print_phylip(os,  valid_indices, line_wrap, true);
			break;
		case eSF_PHYLIP_NONINTERLEAVED:
			print_phylip(os,  valid_indices, line_wrap, false);
			break;
		case eSF_TABDELIMITED:
			{
				// indices range
				cu_iter	index = valid_indices.begin();
				cu_iter	end = valid_indices.end();

				for (; index  != end; ++index)
						os << accession(*index) << "\t" << sequence(*index) << "\n";

			}
			break;
	}

}

}
